{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "gXXhctqjgXO7"
   },
   "source": [
    "##### Copyright 2020 The Cirq Developers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "id": "z2RJVa8qgXou"
   },
   "outputs": [],
   "source": [
    "# @title Licensed under the Apache License, Version 2.0 (the \"License\");\n",
    "# you may not use this file except in compliance with the License.\n",
    "# You may obtain a copy of the License at\n",
    "#\n",
    "# https://www.apache.org/licenses/LICENSE-2.0\n",
    "#\n",
    "# Unless required by applicable law or agreed to in writing, software\n",
    "# distributed under the License is distributed on an \"AS IS\" BASIS,\n",
    "# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.\n",
    "# See the License for the specific language governing permissions and\n",
    "# limitations under the License."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "EQvWLKKRgZR9"
   },
   "source": [
    "# Simulation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "EvZ_JecKga2p"
   },
   "source": [
    "<table class=\"tfo-notebook-buttons\" align=\"left\">\n",
    "  <td>\n",
    "    <a target=\"_blank\" href=\"https://quantumai.google/cirq/simulate/simulation\"><img src=\"https://quantumai.google/site-assets/images/buttons/quantumai_logo_1x.png\" />View on QuantumAI</a>\n",
    "  </td>\n",
    "  <td>\n",
    "    <a target=\"_blank\" href=\"https://colab.research.google.com/github/quantumlib/Cirq/blob/main/docs/simulate/simulation.ipynb\"><img src=\"https://quantumai.google/site-assets/images/buttons/colab_logo_1x.png\" />Run in Google Colab</a>\n",
    "  </td>\n",
    "  <td>\n",
    "    <a target=\"_blank\" href=\"https://github.com/quantumlib/Cirq/blob/main/docs/simulate/simulation.ipynb\"><img src=\"https://quantumai.google/site-assets/images/buttons/github_logo_1x.png\" />View source on GitHub</a>\n",
    "  </td>\n",
    "  <td>\n",
    "    <a href=\"https://storage.googleapis.com/tensorflow_docs/Cirq/docs/simulate/simulation.ipynb\"><img src=\"https://quantumai.google/site-assets/images/buttons/download_icon_1x.png\" />Download notebook</a>\n",
    "  </td>\n",
    "</table>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "bd9529db1c0b"
   },
   "outputs": [],
   "source": [
    "try:\n",
    "    import cirq\n",
    "except ImportError:\n",
    "    print(\"installing cirq...\")\n",
    "    !pip install --quiet cirq\n",
    "    print(\"installed cirq.\")\n",
    "    import cirq"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "4d599188f678"
   },
   "source": [
    "Cirq comes with built-in Python simulators for testing\n",
    "small circuits. The two main types of simulations that Cirq\n",
    "supports are pure state and mixed state. The pure state simulations\n",
    "are supported by `cirq.Simulator` and the mixed state\n",
    "simulators are supported by `cirq.DensityMatrixSimulator`.\n",
    "\n",
    "The names *pure state simulator* and *mixed state\n",
    "simulators* refer to the fact that these simulations are\n",
    "for quantum circuits; including unitary, measurements, and noise\n",
    "that keeps the evolution in a pure state (i.e. a single quantum state)\n",
    "or a mixed state (a mix of quantum states, each with a classical\n",
    "probability). Noisy evolutions are supported by the pure state\n",
    "simulator, as long as they preserve the purity of the state. If you are interested in truly noisy simulation that may not preserve purity, see the [Noisy Simulation](./noisy_simulation.ipynb) page. \n",
    "\n",
    "Some external high-performance simulators also provide an interface\n",
    "to Cirq. These can often provide results faster than Cirq's\n",
    "built-in simulators, especially when working with larger circuits.\n",
    "For details on these tools, see the\n",
    "[external simulators section](#external-simulators).\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "de8b974f6597"
   },
   "source": [
    "## Introduction to pure state simulation\n",
    "\n",
    "Here is a simple circuit:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "f3c82c1c12a6"
   },
   "outputs": [],
   "source": [
    "q0 = cirq.GridQubit(0, 0)\n",
    "q1 = cirq.GridQubit(1, 0)\n",
    "\n",
    "\n",
    "def basic_circuit(meas=True):\n",
    "    sqrt_x = cirq.X**0.5\n",
    "    yield sqrt_x(q0), sqrt_x(q1)\n",
    "    yield cirq.CZ(q0, q1)\n",
    "    yield sqrt_x(q0), sqrt_x(q1)\n",
    "    if meas:\n",
    "        yield cirq.measure(q0, key='q0'), cirq.measure(q1, key='q1')\n",
    "\n",
    "\n",
    "circuit = cirq.Circuit()\n",
    "circuit.append(basic_circuit())\n",
    "\n",
    "print(circuit)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "b4e74a859e3a"
   },
   "source": [
    "You can simulate this by creating a ``cirq.Simulator`` and\n",
    "passing the circuit into its ``run`` method:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "5422e0c06102"
   },
   "outputs": [],
   "source": [
    "simulator = cirq.Simulator()\n",
    "result = simulator.run(circuit)\n",
    "\n",
    "print(result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "31212f98dfa7"
   },
   "source": [
    "The method `run()` returns a\n",
    "[`Result`](https://github.com/quantumlib/Cirq/blob/dd1b1a4a0e26775bf94eb567591397f1f989ca55/cirq-core/cirq/study/result.py#L85).\n",
    "As you can see, the object `result` contains the result of any\n",
    "measurements for the simulation run.\n",
    "\n",
    "The actual measurement results depend on the `random` seed generator (`numpy`).\n",
    "You can set this seed by passing an integer or `numpy.RandomState` as the\n",
    "`seed` parameter in the `Simulator` constructor.\n",
    "\n",
    "Another run can result in a different set of measurement results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "e83b165ebf06"
   },
   "outputs": [],
   "source": [
    "result = simulator.run(circuit)\n",
    "\n",
    "print(result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "b6847671058f"
   },
   "source": [
    "The `run()` methods (`run()` and `run_sweep()`) are designed to mimic what\n",
    "running a program on a quantum computer is actually like. `result` only\n",
    "contains measurement data, and the complete state vector is hidden.\n",
    "\n",
    "### Accessing the state vector\n",
    "\n",
    "To access the full state vector, the `simulate()` methods\n",
    "(`simulate()`, `simulate_sweep()`, `simulate_moment_steps()`)\n",
    "can be used instead. This behavior is only possible in\n",
    "simulation, but can be useful for debugging a circuit:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "e085289e0002"
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "circuit = cirq.Circuit()\n",
    "circuit.append(basic_circuit(False))\n",
    "result = simulator.simulate(circuit, qubit_order=[q0, q1])\n",
    "\n",
    "print(np.around(result.final_state_vector, 3))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "jN69PPL8Tv3L"
   },
   "source": [
    "`simulate()` returns a `SimulationTrialResult` containing\n",
    "the final state, as seen above. The built-in Cirq simulator returns a\n",
    "[`StateVectorTrialResult`](https://github.com/quantumlib/Cirq/blob/dd1b1a4a0e26775bf94eb567591397f1f989ca55/cirq-core/cirq/sim/state_vector_simulator.py#L147)\n",
    ",\n",
    "which includes a number of utilities for analyzing the final state\n",
    "vector.\n",
    "\n",
    "Note that the simulator uses numpy's `float32` precision\n",
    "(which is `complex64` for complex numbers) by default,\n",
    "but that the simulator can take in a  `dtype` of `np.complex128`\n",
    "if higher precision is needed."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "k3hqZBKXVtPK"
   },
   "source": [
    "### Expectation values\n",
    "\n",
    "For applications that measure expectation values of observables,\n",
    "the `simulate_expectation_values()` method provides a simple\n",
    "interface for returning just the desired expectation values.\n",
    "This can be more efficient than returning the entire state vector,\n",
    "particularly when handling multiple results at once, or when using\n",
    "an external simulator."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "mjUGQ78IiQAG"
   },
   "outputs": [],
   "source": [
    "XX_obs = cirq.X(q0) * cirq.X(q1)\n",
    "ZZ_obs = cirq.Z(q0) * cirq.Z(q1)\n",
    "ev_list = simulator.simulate_expectation_values(\n",
    "    cirq.Circuit(basic_circuit(False)), observables=[XX_obs, ZZ_obs]\n",
    ")\n",
    "print(ev_list)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "BL6Ij97Zj7rd"
   },
   "source": [
    "`simulate_expectation_values()` returns a list of\n",
    "expectation values, one for each observable provided."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "028defe596bc"
   },
   "source": [
    "## Qubit and Amplitude Ordering\n",
    "\n",
    "The `qubit_order` argument to the simulator's `run()` method\n",
    "determines the ordering of some results, such as the\n",
    "amplitudes in the final wave function. The `qubit_order` argument is optional: when it is omitted, qubits are ordered\n",
    "ascending by their name (i.e., what `str(qubit)` returns).\n",
    "\n",
    "The simplest `qubit_order` value you can provide is a list of\n",
    "the qubits in the desired order. Any qubits from the circuit\n",
    "that are not in the list will be ordered using the\n",
    "default `str(qubit)` ordering, but come after qubits that are in\n",
    "the list. Be aware that all qubits in the list are included in\n",
    "the simulation, even if they are not operated on by the circuit.\n",
    "\n",
    "The mapping from the order of the qubits to the order of the\n",
    "amplitudes in the wave function can be tricky to understand.\n",
    "Basically, it is the same as the ordering used by `numpy.kron`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "982a827446d3"
   },
   "outputs": [],
   "source": [
    "outside = [1, 10]\n",
    "inside = [1, 2]\n",
    "print(np.kron(outside, inside))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "a36137200917"
   },
   "source": [
    "More concretely, the `k`'th amplitude in the wave function\n",
    "will correspond to the `k`'th case that would be encountered\n",
    "when nesting loops over the possible values of each qubit.\n",
    "\n",
    "The first qubit's computational basis values are looped over\n",
    "in the outermost loop, the last qubit's computational basis\n",
    "values are looped over in the inner-most loop, etc.:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "8141e8461313"
   },
   "outputs": [],
   "source": [
    "i = 0\n",
    "for first in [0, 1]:\n",
    "    for second in [0, 1]:\n",
    "        print('amps[{}] is for first={}, second={}'.format(i, first, second))\n",
    "        i += 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "48c4d1acd296"
   },
   "source": [
    "You can check that this is in fact the ordering with a\n",
    "circuit that flips one qubit out of two:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "befc299bfd2a"
   },
   "outputs": [],
   "source": [
    "q_stay = cirq.NamedQubit('q_stay')\n",
    "q_flip = cirq.NamedQubit('q_flip')\n",
    "c = cirq.Circuit(cirq.X(q_flip))\n",
    "\n",
    "# first qubit in order flipped\n",
    "result = simulator.simulate(c, qubit_order=[q_flip, q_stay])\n",
    "print(abs(result.final_state_vector).round(3))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "ca70c50a782f"
   },
   "outputs": [],
   "source": [
    "# second qubit in order flipped\n",
    "result = simulator.simulate(c, qubit_order=[q_stay, q_flip])\n",
    "print(abs(result.final_state_vector).round(3))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "deda90b07995"
   },
   "source": [
    "## Stepping through a circuit\n",
    "\n",
    "When debugging, it is useful to not just see the end\n",
    "result of a circuit, but to inspect the state of the system\n",
    "at different steps in the circuit.  \n",
    "\n",
    "To support this, Cirq provides a method to return an iterator\n",
    "over a `Moment` by `Moment` simulation.  This method is named `simulate_moment_steps`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "68bc0a762450"
   },
   "outputs": [],
   "source": [
    "circuit = cirq.Circuit()\n",
    "circuit.append(basic_circuit())\n",
    "for i, step in enumerate(simulator.simulate_moment_steps(circuit)):\n",
    "    print('state at step %d: %s' % (i, np.around(step.state_vector(copy=True), 3)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "ef682175a07b"
   },
   "source": [
    "The object returned by the `moment_steps` iterator is a\n",
    "`StepResult`. This object has the state along with any\n",
    "measurements that occurred before or during that step."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "f4298ecad1ef"
   },
   "source": [
    "### Alternate stepping behavior\n",
    "\n",
    "For simulators that do not support `simulate_moment_steps`,\n",
    "it is possible to replicate this behavior by splitting\n",
    "the circuit into \"chunks\" and passing the results of each\n",
    "chunk as the initial state for the next chunk:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "8c8c037355eb"
   },
   "outputs": [],
   "source": [
    "chunks = [cirq.Circuit(moment) for moment in basic_circuit()]\n",
    "next_state = 0  # represents the all-zero state\n",
    "for i, chunk in enumerate(chunks):\n",
    "    result = simulator.simulate(chunk, initial_state=next_state)\n",
    "    next_state = result.final_state_vector\n",
    "    print(f'state at step {i}: {np.around(next_state, 3)}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "ca2bc7e8b96b"
   },
   "source": [
    "The added cost of passing state vectors around like this\n",
    "is nontrivial; for this reason, this workaround should\n",
    "only be used with simulators that do not support\n",
    "`simulate_moment_steps`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "p_IOhhsWgnTg"
   },
   "source": [
    "## Parameterized values and studies\n",
    "\n",
    "In addition to circuit gates with fixed values, Cirq also\n",
    "supports gates which can have `Symbol` values (see\n",
    "[Gates](../build/gates.ipynb)). These are values that can be resolved\n",
    "at runtime. \n",
    "\n",
    "For simulators, these values are resolved by\n",
    "providing a `cirq.ParamResolver`.  A `cirq.ParamResolver` provides\n",
    "a map from the `Symbol`'s name to its assigned value."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "1b5c6f4438ef"
   },
   "outputs": [],
   "source": [
    "import sympy\n",
    "\n",
    "rot_w_gate = cirq.X ** sympy.Symbol('x')\n",
    "circuit = cirq.Circuit()\n",
    "circuit.append([rot_w_gate(q0), rot_w_gate(q1)])\n",
    "print(circuit)\n",
    "for y in range(5):\n",
    "    resolver = cirq.ParamResolver({'x': y / 4.0})\n",
    "    result = simulator.simulate(circuit, resolver)\n",
    "    print(f\"params:{result.params}, state vector:{np.round(result.final_state_vector, 2)}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "10617793bfd2"
   },
   "source": [
    "In the previous example, the symbol `x` is used in two gates, and then the resolver\n",
    "provides this value at run time.\n",
    "\n",
    "Parameterized values are most useful in defining what is called a\n",
    "\"sweep\", which is a sequence of trials, where each\n",
    "trial is a run with a particular set of parameter values.\n",
    "\n",
    "Running a sweep returns a `Result` for each set of fixed parameter\n",
    "values and repetitions.  \n",
    "\n",
    "For instance:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "ccc88952ba49"
   },
   "outputs": [],
   "source": [
    "resolvers = [cirq.ParamResolver({'x': y / 2.0}) for y in range(3)]\n",
    "circuit = cirq.Circuit()\n",
    "circuit.append([rot_w_gate(q0), rot_w_gate(q1)])\n",
    "circuit.append([cirq.measure(q0, key='q0'), cirq.measure(q1, key='q1')])\n",
    "results = simulator.run_sweep(program=circuit, params=resolvers, repetitions=2)\n",
    "for result in results:\n",
    "    print(f\"params:{result.params}\")\n",
    "    print(f\"measurements:\")\n",
    "    print(result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "5ca5f417ca71"
   },
   "source": [
    "The previous example demonstrates that assigning different values to gate parameters yields\n",
    "different results for each trial in the sweep, and that each trial is repeated\n",
    "`repetitions` times.  See [Parameter Sweeps](./params.ipynb) for more information on sweeping and parameters."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "d449810f2cad"
   },
   "source": [
    "## Mixed state simulations\n",
    "\n",
    "In addition to pure state simulation, Cirq also supports\n",
    "simulation of mixed states. \n",
    "\n",
    "Even though this simulator is not as efficient as the pure state simulators,\n",
    "they allow for a larger class of noisy circuits to be run as well as keeping\n",
    "track of the simulation's density matrix. This fact can allow for more exact\n",
    "simulations: the density matrix can represent all possible results of a\n",
    "noisy circuit, while the pure-state simulator can only sample from these\n",
    "results.\n",
    "\n",
    "Mixed state simulation is supported by the\n",
    "`cirq.DensityMatrixSimulator` class.\n",
    "\n",
    "Here is a simple example of simulating a channel using the\n",
    "mixed state simulator:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "9bc8775b8ef9"
   },
   "outputs": [],
   "source": [
    "q = cirq.NamedQubit('a')\n",
    "circuit = cirq.Circuit(cirq.H(q), cirq.amplitude_damp(0.2)(q), cirq.measure(q))\n",
    "simulator = cirq.DensityMatrixSimulator()\n",
    "result = simulator.run(circuit, repetitions=100)\n",
    "print(result.histogram(key='a'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "164ccf2b7d71"
   },
   "source": [
    "The previous example creates a state in an equal superposition of 0 and 1, then applies amplitude damping which takes 1 to 0 with something like a probability of 0.2. \n",
    "\n",
    "You can see that, instead of about 50 percent of the timing being in 0, about 20 percent of the 1 has been converted into 0, so you end up with total around 60 percent in the 0 state.\n",
    "\n",
    "Like the pure state simulators, the mixed state simulator supports `run()` and `run_sweeps()` methods. \n",
    "\n",
    "The `cirq.DensityMatrixSimulator` also supports getting access to the density matrix of the circuit at the end of simulating the circuit, or when stepping through the circuit.  These are done by the `simulate()` and `simulate_sweep()` methods, or, for stepping through the circuit, via the `simulate_moment_steps` method. For example, you can simulate creating an equal superposition followed by an amplitude damping channel with a gamma of 0.2 by:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "c66bc5094882"
   },
   "outputs": [],
   "source": [
    "q = cirq.NamedQubit('a')\n",
    "circuit = cirq.Circuit(cirq.H(q), cirq.amplitude_damp(0.2)(q))\n",
    "simulator = cirq.DensityMatrixSimulator()\n",
    "result = simulator.simulate(circuit)\n",
    "print(np.around(result.final_density_matrix, 3))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "ce9b181d8c72"
   },
   "source": [
    "Note: The density matrix at the end of the simulation is accessible via `final_density_matrix`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "3e15de52ca5f"
   },
   "source": [
    "## External simulators\n",
    "\n",
    "There are a few high-performance circuit simulators which\n",
    "provide an interface for simulating Cirq `Circuit`s.\n",
    "These projects are listed below, along with their PyPI package\n",
    "name and a description of simulator methods that they support.\n",
    "\n",
    "For most users we recommend [qsim](https://github.com/quantumlib/qsim),\n",
    "but each simulator is optimized for specific use cases. Before choosing\n",
    "a simulator, make sure it supports the behavior that you need!\n",
    "\n",
    "| Project name | PyPI package | Description |\n",
    "| --- | --- | --- |\n",
    "| [qsim](https://github.com/quantumlib/qsim) | qsimcirq | Implements `cirq.SimulatesAmplitudes`, `cirq.SimulatesFinalState`, and `cirq.SimulatesExpectationValues`. Recommended for deep circuits with up to 30 qubits (consumes 8GB RAM). Larger circuits are possible, but RAM usage doubles with each additional qubit. |\n",
    "| [qsimh](https://github.com/quantumlib/qsim/blob/main/qsimcirq/qsimh_simulator.py) | qsimcirq | Implements `cirq.SimulatesAmplitudes`. Intended for heavy parallelization across several computers; Cirq users should generally prefer qsim. |\n",
    "| [quimb](https://quimb.readthedocs.io/en/latest/) | quimb | Not based on `cirq.Simulator`; instead uses a Cirq-to-quimb translation layer provided in `contrib/quimb`. In addition to circuit simulation, this allows the use of quimb circuit-analysis tools on Cirq circuits. |\n",
    "| [qFlex](https://github.com/ngnrsaa/qflex) | qflexcirq | Implements `cirq.SimulatesAmplitudes`. RAM usage is highly dependent on the number of two-qubit gates in the circuit. Not recommended - prefer qsim or quimb. |\n"
   ]
  }
 ],
 "metadata": {
  "colab": {
   "name": "simulation.ipynb",
   "toc_visible": true
  },
  "kernelspec": {
   "display_name": "Python 3",
   "name": "python3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
